#ifndef POLYNOMIAL_H
#define POLYNOMIAL_H

#include <iostream>
#include <vector>
#include <cmath>
#include <algorithm>

using namespace std;

// Declare class
template <typename T>
class Polynomial;
// Declare friend functions
template <typename T>
ostream& operator<<(ostream& os, const Polynomial<T>& p);
template <typename T>
const Polynomial<T> operator+(const Polynomial<T>& p1, const Polynomial<T>& p2);
template <typename T>
const Polynomial<T> operator-(const Polynomial<T>& p1, const Polynomial<T>& p2);
template <typename T>
const Polynomial<T> operator*(const Polynomial<T>& p1, const Polynomial<T>& p2);
template <typename T>
const Polynomial<T> operator*(const T x, const Polynomial<T>& p);

template <typename T>
class Polynomial
{
private:
    vector<T> coefficients;
public:
    Polynomial();
    Polynomial(const vector<T> vec);
    ~Polynomial();

    // Set coefficients of polynomial
    void setCoe(const vector<T> vec);
    // Adjust coefficients and ensure leading coefficient not zero
    void adjustCoe();
    // Get coefficients of polynomial
    const vector<T>& getCoe() const;
    // Get degree of polynomial
    const int getDegree() const;
    // Get value of polynomial at x
    const T getValue(const T& x) const;
    // Get derivation of Polynomial
    Polynomial<T> getDerivation() const;
    
    // Overload operators used in polynomials
    friend ostream& operator<< <> (ostream& os, const Polynomial<T>& p);
    friend const Polynomial<T> operator+ <> (const Polynomial<T>& p1, const Polynomial<T>& p2);
    friend const Polynomial<T> operator- <> (const Polynomial<T>& p1, const Polynomial<T>& p2);
    friend const Polynomial<T> operator* <> (const Polynomial<T>& p1, const Polynomial<T>& p2);
    friend const Polynomial<T> operator* <> (const T x, const Polynomial<T>& p);
};

template <typename T>
Polynomial<T>::Polynomial()
{
    double d[1] = {0};
    vector<T> coe(d,d+1);
    coefficients = coe;
}

template <typename T>
Polynomial<T>::Polynomial(const vector<T> vec)
{
    setCoe(vec);
    adjustCoe();
}

template <typename T>
Polynomial<T>::~Polynomial()
{};

template <typename T>
void Polynomial<T>::setCoe(const vector<T> vec)
{
    coefficients = vec;
}

template <typename T>
void Polynomial<T>::adjustCoe()
{
    while(coefficients.back() == 0 && coefficients.size() > 1)
	coefficients.pop_back();
}

template <typename T>
const vector<T>& Polynomial<T>::getCoe() const
{
    return coefficients;
}

template <typename T>
const int Polynomial<T>::getDegree() const
{
    if (coefficients.size() > 1)
	return coefficients.size() - 1;
    else if (coefficients[0] == 0)
	return -1;
    else
	return 0;
}

template <typename T>
const T Polynomial<T>::getValue(const T& x) const
{
    T sum = 0;
    const int degree = this->getDegree();
    for (int i = 0; i <= degree; ++i)
	sum += coefficients[i] * pow(x,i);
    return sum;
}

template <typename T>
Polynomial<T> Polynomial<T>::getDerivation() const
{
    const vector<T>& vec = this->getCoe();
    const int degree = this->getDegree();
    if (degree > 0)
    {
	vector<T> dvec(degree);
	for (int i = 0; i < degree; ++i)
	    dvec[i] = vec[i+1] * (i+1);
	Polynomial<T> l1(dvec);
	return l1;
    }
    else
    {
	Polynomial<T> l;
	return l;
    }
}

template <typename T>
ostream& operator<<(ostream& os, const Polynomial<T>& p)
{
  vector<double> v = p.getCoe();   
  const int degree = p.getDegree();
  if (degree == -1){
    os << 0;
    return os;
  }
  else if (degree == 0){
       os << v[0];
       return os;
  }
  else{
     if (v[degree] > 0 && v[degree] != 1)
          os << v[degree] << "x^" << degree;
     else if (v[degree] < 0 && v[degree] != -1)
          os << v[degree] << "x^" << degree;
     else if (v[degree] == 1)
          os << "x^" << degree;
     else if (v[degree] == -1)
          os << "-x^" << degree;
    for (int i = degree - 1 ; i > 0 ; i--){
      if (v[i] > 0 && v[i] != 1)
           os << "+" << v[i] << "x^" << i;
      else if (v[i] < 0 && v[i] != -1)
           os << v[i] << "x^" << i;
      else if (v[i] == 1)
           os << "+" << "x^" << i;
      else if (v[i] == -1)
           os << "-x^" << i;
    }
    if (v[0] > 0)
         os << "+" << v[0];
    else if (v[0] < 0)
         os << v[0];
    return os;
  }
    
}

template <typename T>
const Polynomial<T> operator+(const Polynomial<T>& p1, const Polynomial<T>& p2)
{
     vector<T> coe1 = p1.getCoe();
     vector<T> coe2 = p2.getCoe();
     const int m = coe1.size();
     const int n = coe2.size();
     const int M = max(m,n);
     vector<T> coe(M);
     coe1.resize(M,0);
     coe2.resize(M,0);
     for (int i = 0 ; i < M ; i++)
          coe[i] = coe1[i] + coe2[i];
     Polynomial<T> p(coe);
     p.adjustCoe();
     return p;
}

template <typename T>
const Polynomial<T> operator-(const Polynomial<T>& p1, const Polynomial<T>& p2)
{
     vector<T> coe1 = p1.getCoe();
     vector<T> coe2 = p2.getCoe();
     const int m = coe1.size();
     const int n = coe2.size();
     const int M = max(m,n);
     vector<T> coe(M);
     coe1.resize(M,0);
     coe2.resize(M,0);
     for (int i = 0 ; i < M ; i++)
          coe[i] = coe1[i] - coe2[i];
     Polynomial<T> p(coe);
     p.adjustCoe();
     return p;
}

template <typename T>
const Polynomial<T> operator*(const Polynomial<T>& p1, const Polynomial<T>& p2)
{
     vector<T> coe1 = p1.getCoe();
     vector<T> coe2 = p2.getCoe();
     const int m = coe1.size();
     const int n = coe2.size();
     const int M = m + n;
     vector<T> coe(M,0);
     coe1.resize(M,0);
     coe2.resize(M,0);
     for (int i = 0 ; i < M ; i++){
          for (int j = 0 ; j <= i ; j ++){
               coe[i] += coe1[j] * coe2[i - j];
          }
     }  
     Polynomial<T> p(coe);
     p.adjustCoe();
     return p;
}

template <typename T>
const Polynomial<T> operator*(const T x, const Polynomial<T>& p)
{
     vector<T> coe = p.getCoe();
     const int m = coe.size();
     vector<T> coen(m,0);
     for (int i = 0 ; i < m ; i++)
          coen[i] = x*coe[i];
     Polynomial<T> p1(coen);
     p1.adjustCoe();
     return p1;
}


#else
//do nothing
#endif
